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We consider the dynamics of a dark soliton in an elongated harmonically 
trapped Bose-Einstein condensate. A central question concerns the behavior 
at finite temperatures, where dissipation arises due to the presence of a ther- 
mal cloud. We study this problem using coupled Gross-Pitaevskii and N-body 
simulations, which include the mean field coupling between the condensate 
and thermal cloud. We find that the soliton decays relatively quickly even at 
very low temperatures, with the decay rate increasing with rising temperature. 
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1. INTRODUCTION 

Dark solitons represent a ubiquitous feature of nonlinear system^, and 
their experimental realization in atomic Bose-Einstein condensates (BECs ^^ 4 * 5 * 
has opened up new possibilities for their study. The stability of the soliton, 
and its dependence upon the geometry of the condensate, is a problem of 
considerable interest. In three dimensions a soliton is dynamically unstable 
to undulations in its transverse profile. This so-called "snake" instability 
leads to the soliton decaying to vortex lines or rings, and has been observed 
experimentally 4 5 . This instability, however, can be suppressed by confining 
the condensate much more strongly in the transverse direction than longitu- 
dinally, so that the condensate is highly elongated along the axial direction. 

A second decay mechanism arises from the fact that at non-zero tem- 
peratures the condensate coexists with a noncondensed cloud composed of 
thermally excited quasiparticles. Interactions of the soliton with the thermal 
cloud leads to dissipation, so that the soliton loses energy and accelerates. 
Experiment^ suggest that the timescale for this decay can be relatively 
short, on the order of 10 ms, so has important consequences for the dynam- 
ics. This decay process has previously been studied theoretically by Fedichev 
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et alW and Muryshev et a/Pi by considering the reflection of thermal excita- 
tions from the soliton. The role of quantum fluctuations was addressed in 
Dziarmaga et 

In this paper we describe the use of numerical simulations 9 to model a 
soliton at finite temperatures in a highly elongated condensate, where the 
trap frequency in the radial direction, u±, greatly exceeds that along the 
axial direction, u> z . At T = the soliton is predicte d 10 * 11 *, on the basis of 
the Gross-Pitaevskii (GP) equation, to oscillate longitudinally at a frequency 
of lu z /^/2, which involves continuous undamped soliton-sound interactional. 
We follow the motion of the soliton as a function of time, finding the expected 
oscillation along the axial direction. At finite temperature the energy loss of 
the soliton due to the presence of the thermal cloud leads to an increase in 
the amplitude of the oscillation. We investigate this soliton decay process as 
a function of temperature, and find that it is appreciable even at very low 
T. 



2. THEORY 

Our treatment of finite temperatures is based on the ZNG formalism^, 
where we solve the following coupled equations 

<9^ / % 2 V 2 \ 
ih —={-—- + V + 9 n c + 2 9 nj% (1) 

^ + --V/-VC/-V p / = 0. (2) 

at m 

The former is a generalized GP equation for the condensate wavefunction 
ty(r,t), and the latter is a Boltzmann equation for the thermal cloud phase 
space density f(p,r,t). Here the condensate and thermal cloud densities 
are defined as n c = \^\ 2 and n = J f dp/h 3 respectively, while g = A7rh 2 a/m 
parametrizes the mean field interactions between atoms of mass m with 
scattering length a. The effective potential is U = V + 2g(n c + n), with V = 
m{id\r 2 -\- uj 2 z 2 ) /2 representing the harmonic trap confining the atoms. Note 
that for the purposes of this paper we neglect collisions between the atoms, 
and the coupling between the condensate and thermal cloud is purely mean 
field in nature. For the parameters considered here, inclusion of collisions is 
expected to slightly increase the damping rate, but this will be considered 
explicitly in future work. 

The initial condition for our simulations is found by self-consistently 
solving in imaginary time to find the ground state wavefunction, which 
is coupled to a thermal cloud described by a Bose distribution f(p,r) = 
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{exp[(p 2 /2m + U - Li)/k B T] - l}" 1 , where /x is the condensate chemical 
potential. The methods for solving Eqs. Q and © for the subsequent finite 
temperature dynamics are described by Jackson and ZarembcP. A difference 
with respect to this earlier work is that here the GP equation (JJ) is solved 
using cylindrical coordinates with a Crank-Nicholson time-stepping scheme. 
Since the system is cylindrically symmetric the problem effectively reduces 
to two dimensions, improving the efficiency of the simulations. 



3. RESULTS 

Finite temperature simulations have been performed for N = 2 x 10 4 
87 Rb atoms in a trap with frequencies u>± = 2ir x 2500 Hz and u> z = 2-ir x 
10 Hz, resulting in a cloud that is highly elongated along the axial direc- 
tion. The critical temperature for Bose condensation is close to that for 
an ideal gas in the thermodynamic limiPl, T c ° = 0.941ft(c4i^A0 1/3 AB, 
which for these parameters corresponds to = 486 nK. At time t = 
the equilibrium configuration for a given temperature is imprinted with 
a dark soliton, by multiplying the condensate wavefunction Vl/(r) witrP^ 
^ s (z) = /3tanh(/3z/£) +iv/c, where = i/l — (v/c)' 2 , £ = h/^/mgn is 
the condensate healing length, and c = y / gn/(2m) is the speed of sound for 
density n. For the benefit of the following discussion, it is useful to note 
that for a condensate of uniform density no the soliton energy is given by^ 
E s = 4ftcno/3 3 /3, while the density at the soliton minimum (equal to the 
difference between the background density and the depth of the soliton) is^ 
n s = tiqv 2 j(? . So an increase in v corresponds to a decrease in the soliton 
energy and depth, which both tend to zero as v — * c. 

For T = the dynamics are simply given by the GP equation Q with 
h = 0. The black lines in Fig. ^pfot the subsequent evolution of the density 
at early times, where the soliton oscillates along the axial direction. The 
position of the soliton minimum, z s , is plotted in Fig. [2(a), illustrating that 
the soliton oscillates at constant amplitude with a frequency of u> z /y2, as 
expected 10 * 11 * 13 *. The parameter n s , representing the density at the soliton 
minimum as labelled in Fig.^( a )) acquires its maximum value at the center 
of the trap and is zero at the turning points of the soliton motion. Fig. |2] (b) 
plots n s , divided by the central density of the condensate without a soliton, 
no- The peak value of this parameter (attained at the trap centre) remains 
constant, as expected for undamped motion. 

The corresponding density cross-sections for T = 150 nK are plotted 
with gray lines in Fig.Q where one sees that, even during the first oscillation, 
the presence of the thermal cloud leads to an increase in the amplitude of 
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Fig. 1. Density |\t| 2 along an axial cross-section of the condensate, showing 
the evolution of a dark soliton with initial velocity v = 0.5c. The black lines 
plot the densities for T = 0, while the gray lines are for T = 150 nK, at times 
of (a) uj z t = 0.2, (b) uo z t = 2.2, (c) u z t = 4.4, and (d) u z t = 6.7. The axial 
position z is in units of a z = y 'fi / '(muj z ) , while the density is in arbitrary 
units. The density at the soliton minimum, n s , is labelled in (a). 

the oscillation (Fig. ^ (b) and (d)), in addition to a significant decrease in 
the depth of the soliton when it is at the trap center (Fig.n(c)). The steady 
rise in oscillation amplitude is also shown in Fig. |2] (a) , while the decrease in 
soliton depth at z s = (corresponding to an increase in n s ) is illustrated in 
Fig. [21(b). This accompanies an increase in velocity at z s = 0, and reflects 
an energy decay as the soliton interacts with the thermal cloud. 

Simulations have also been performed at other temperatures, and one 
finds less damping at lower T (as illustrated for T = 100 nK in Fig. |2j) 
and more damping at higher T, as expected. However, it is remarkable 
that appreciable damping occurs even for these relatively low tempera- 
tures (T < 0.3T C ) , which contrasts with the situation for collective modes 
where the effects of finite temperature tend to be significant only at higher 
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Fig. 2. (a) Axial position (in units of a z ) of the soliton minimum as a 
function of time, where initially v = 0.5c. The curves show results for 
different temperatures, with T = (solid black), T = 100 nK (dashed), and 
T = 150 nK (gray), (b) Density at the soliton minimum n s (divided by the 
central density of the condensate without a soliton, no) as a function of time, 
where the curves are labelled as in (a). 

temperatures^EHini Another noteworthy feature relates to the time de- 
pendence of the amplitude, the increase of which appears to saturate as it 
approaches the radius of the condensate. This effect has also been observed 
for v = 0.25c, where the initial oscillation amplitude is smaller, and hence 
the saturation does not occur until later times. 



4. CONCLUSIONS 

We have performed numerical simulations to study the dynamics of 
dark solitons in a Bose-Einstein condensate at different temperatures. Our 
numerics for zero temperature reproduce soliton oscillations along the axial 
direction at constant amplitude and at a frequency of u> z /y2. However, at 
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finite temperature, mean field coupling between the condensate and ther- 
mal cloud is shown to lead to a lower frequency oscillation with a steadily 
increasing amplitude. This corresponds to the soliton losing energy to the 
thermal cloud, with a rate that increases with rising initial temperature. 

Note that this energy transfer cannot be considered as a heating of the 
thermal cloud, since temperature is not defined during the simulation due 
to a lack of rethermalizing collisions. Future calculations will study the in- 
fluence of these collisional processes, as well as modelling the soliton decay 
for the experimental parameters of Burger et alP. A further motivation for 
future work is to ascertain the feasibility of Proukakis et al™, where a para- 
metric driving scheme was proposed for stabilization of a soliton against a 
phenomenologically modelled thermal cloud. We plan to revisit this problem 
with our more detailed finite temperature simulations. 
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